The present invention, in some embodiments thereof, relates to image reconstruction of image data from a plurality of radiation sources and, more particularly, but not exclusively, to image reconstruction for computed tomography (CT).
Methods for imaging an examination object with a CT system are generally known. Projections emitted by the X-ray source can be shaped into a two dimensional (2D) fan shaped beam or a three dimensional (3D) cone shaped beam. Radiation emitted by the X-ray source is attenuated by the examination object and detected by the detector. In known methods an X-ray source and an opposing detector rotate around the examination object while the detector acquires data from a series of projections captured at different angles with respect to the examination object. Output from the detector is processed by reconstruction methods to provide sectional images through the examination object.
Filter Back-Projection (FBP) algorithm is a known reconstruction method used for reconstructing fan beam acquired data. The reconstruction algorithm used for reconstructing cone beam acquired data is usually of the type called “Feldkamp” or “FDK” method or modifications of the FDK method. The FDK method was published by L. A. Feldkamp, L. C. Davis, and J. W. Kress, JOSA A, Vol. 1, Issue 6, pp. 612-619 (1984), the content of which incorporated herein by reference. The reconstructed images obtained by these methods are known to have artifacts that worsen for cone beam angles that are distanced from the mid plane (the source trajectory plane) due to “data incompleteness”. For wide angle cone beams, the artifacts that are produced are known to be significant. Similar problems occur with other cone beam reconstruction algorithms known in the art.
It is also known to scan with a plurality of axial offset X-ray sources or to perform a plurality of axially offset scans with a single source, both for increasing the completeness of the data acquired. Typically, the axial offset provides overlap between the projected cone beams. For overlapping cone beams, a weighting method is required to combine the data obtained from the difference sources. Weighting can be performed on the data captured from the detector on a per view basis during image reconstruction or after image reconstruction, on data acquired from each of the sources. An exemplary dual source CT scanner and method is described for example in U.S. Pat. No. 7,869,561 entitled “CONE-BEAM CT,” assigned to Arineta Ltd., the content of which is incorporated herein by reference.
A publication entitled “A new method to combine 3D reconstruction volumes for multiple parallel circular cone beam orbits” by J. Baek and N. J. Pelc that was published in Medical physics, vol. 37, no. 10 pp. 5351-5360 (2010) the content of which is incorporated herein by reference, describes a method for reconstructing data from a plurality of X-ray sources. It is described that the cone beam projection data of each orbit are separately reconstructed using the FDK algorithm. Subsequently, overlap regions of the reconstructed image volumes are combined using weighted averaging in frequency space. A smoothly varying weighting function in the overlap region is used to avoid image artifacts caused by the abrupt transition in frequency space.
One known disadvantage of applying weighting after image reconstruction is that only volume elements (voxels) that have sufficient angular coverage by a given source can provide an image voxel for that source. For example, it is known that FBP reconstruction can be applied only to voxels that are irradiated from at least 180° plus the beam fan angle. As a result, some of the acquired data is not used in reconstruction and dose efficiency of the scanner is reduced. FIG. 6C shows an exemplary image reconstructed by weighting an image formed separately from each source. Reconstruction by this method can also be seen to include artifacts 651.
U.S. Pat. No. 6,996,204 entitled “Sequential computed tomography method,” the contents of which is incorporated herein by reference describes a computed tomography method in which an examination zone is irradiated over a full rotation of the source and detector, from two mutually offset, preferably circular source trajectories. In an intermediate region the absorption distribution is reconstructed by combining measurement values from both trajectories with weights. The weight for a particular voxel is a function of a distance between the voxel to be reconstructed and the relevant trajectory. The weight is increased as the distance between the voxel to be reconstructed and the relevant trajectory decreases. This method requires obtaining scan data over at least 360 degrees.
A publication entitled “3D analytic cone-beam reconstruction for multiaxial CT acquisitions” by Z. Yin, B. De Man, and J. Pack that was published in International journal of Biomedical imaging, vol. 2009, ID 538389, (2009), the content of which is incorporated herein by reference, describes a methods for reconstruction data from multiple axial sources for both full 360 degrees and partial scans 180 degrees plus fan angle. For the partial scan geometry, the described method uses a view based weighting in order to combine the data from the different sources during reconstruction. The weighting described is determined according to the cone angle between each voxel and each source, with feathering at the borders between regions in order to avoid weight discontinuities.
FIG. 2 shows a simplified schematic illustration of a cross sectional view of two radiation sources irradiating a volume of interest (VOI) and a corresponding example of a prior art linear weighting function used for image reconstructing. In some known methods, the weight of data obtained from each source S1 and S2 at each voxel in a VOI 150 is calculated according to its cone angle or distance from a source plane of S1 and S2. Typically, VOI 150 includes a region A irradiated only by S1, a region B irradiated only by S2 and a region C irradiated by both S1 and S2. Typically, regions A, B and C change with view angle of the scanner. In some known methods, the weighting function for each of S1 and S2 within region C is defined as a linearly varying function that depends on cone angle or distance from the source plane while the weighting function in each of regions A and B is required to be equal to 0 or to 1. Typically, the slope is defined such that the weighting functions at source planes 101 and 102 equal 1 or 0. Since, source planes 101 and 102 are typically within regions A and B respectively and displaced from the boundaries, e.g. ZAC and ZBC that separate regions A, B and C, there are discontinuities in the weighting function at the boundaries and/or the crossover between region C and each of regions A and B. The sudden transitions in the weighting function are known to introduce image artifacts. It is known to use feathering to reduce these image artifacts.
In the example shown in FIG. 2, the weights associated with each source are plotted along the Z direction at two different Y voxel coordinates. The weighting functions 161,162 and 171,172 represent weighting functions calculated as a function of cone angle or distance from the source planes 101, 102 as is known in the art. Discontinuities in the weighting functions at the boundary crossing are shown. Typically, these transitions are more pronounced as the distance from detector 106 increases. For example, the transitions in weighting functions 161 and 162 are more pronounced than the transitions in weighting functions 171 and 172 defined for a plane that is closer to detector 106. Known methods apply feathering, e.g. feathering 161′, 162′, 171′, 172′. FIG. 2 also shows exemplary weighting functions 181, 182 associated with each source S1, S2 at an exemplary Z coordinate and along a Y direction. Typically, for linearly varying weight functions along the Z direction, the weight functions along the Y direction are constant. Typically, feathering is also applied along the Y direction.
An example of image artifacts that develop due to discontinuity of the weighting function at boundary regions, e.g. when feathering is not applied is shown in FIG. 6A. Typically, when feathering is not applied, the discontinuity of the weighting function leads to artifacts, e.g. artifact 650 in FIG. 6A. FIG. 6B shows an exemplary image reconstructed with feathering. Typically, artifact 650 is less pronounced although still apparent when feathering is applied.